package br.edu.ufcg.msnlab2.grupo02.booleSimpson;

import br.edu.ufcg.msnlab2.misc.Function;
import br.edu.ufcg.msnlab2.misc.NumericalIntegrationMethod;

public class SimpsonsAdaptativeRule implements NumericalIntegrationMethod {

	private static final double DEFAULT_TOLERANCE = 0.01;
	private static final int DEFAULT_MAXIMUM_RECURSION_DEPTH = 40;
	
	private final double tolerance;
	private final int maximumRecursionDepth;

	public SimpsonsAdaptativeRule() {
		this(DEFAULT_TOLERANCE);
	}
	
	public SimpsonsAdaptativeRule( double tolerance ) {
		this( tolerance, DEFAULT_MAXIMUM_RECURSION_DEPTH );
	}
	
	public SimpsonsAdaptativeRule( double tolerance, int maximumRecursionDepth ) {
		this.tolerance = tolerance;
		this.maximumRecursionDepth = maximumRecursionDepth;
	}

	public double adaptiveSimpsonsAux( Function f, double a, double b, double epsilon, double S, double fa, double fb, double fc, int bottom ) {
		double c = ( a + b ) / 2, h = b - a;
		double d = ( a + c ) / 2, e = ( c + b ) / 2;
		double fd = f.evaluate( d );
		double fe = f.evaluate( e );
		double Sleft = ( h / 12 ) * ( fa + 4 * fd + fc );
		double Sright = ( h / 12 ) * ( fc + 4 * fe + fb );
		double S2 = Sleft + Sright;
		if ( bottom <= 0 || f.evaluate( Math.abs( S2 - S ) ) <= 15 * epsilon )
			return S2 + ( S2 - S ) / 15;
		return adaptiveSimpsonsAux( f, a, c, epsilon / 2, Sleft, fa, fc, fd, bottom - 1 ) + adaptiveSimpsonsAux( f, c, b, epsilon / 2, Sright, fc, fb, fe, bottom - 1 );
	}

	@Override
	public double integrate( Function f, double a, double b ) {
		
		double c = ( a + b ) / 2;
		double h = b - a;
		double fa = f.evaluate( a );
		double fb = f.evaluate( b );
		double fc = f.evaluate( c );
		double S = ( h / 6 ) * ( fa + 4 * fc + fb );
		return adaptiveSimpsonsAux( f, a, b, tolerance, S, fa, fb, fc, maximumRecursionDepth );
	}

}
